To quantitatively examine the efficacy of vegetation restoration in drylands globally.
Study-level viz to document patterns in exclusions primarily and the relatie frequenices, at the study level, of major categories of evidence.
#study data####
library(tidyverse)
studies <- read_csv("data/studies.csv")
studies
## # A tibble: 278 x 18
## ID title technique data region exclude rationale observations
## <dbl> <chr> <chr> <chr> <chr> <chr> <chr> <chr>
## 1 152 Shor… seeding,… expe… Africa no <NA> <NA>
## 2 180 Rest… chemical… App.… Africa no <NA> <NA>
## 3 229 Infl… soil see… fiel… Africa no <NA> <NA>
## 4 230 Acti… planting fiel… Africa no <NA> <NA>
## 5 255 The … grazing … fiel… Africa no <NA> <NA>
## 6 262 Reve… seeding,… eper… Africa no <NA> <NA>
## 7 263 The … phytogen… fiel… Africa no <NA> <NA>
## 8 264 Eval… seeding,… fiel… Africa no <NA> <NA>
## 9 271 Patc… natural … fiel… Africa no <NA> <NA>
## 10 4 Fact… natural … App.… Africa no <NA> <NA>
## # ... with 268 more rows, and 10 more variables: disturbance <chr>,
## # system <chr>, goal <chr>, intervention <chr>, paradigm <chr>,
## # grazing <chr>, hypothesis <chr>, soil <chr>, benchmark <chr>,
## # notes <chr>
#quick look at rationale needed
exclusions <- studies %>%
filter(exclude == "yes")
#quick look at studies with paradigms
evidence <- studies %>%
filter(exclude == "no")
#library(skimr)
#skim(evidence)
#study-level viz#####
#exclusions
#ggplot(exclusions, aes(rationale, fill = region)) +
# geom_bar() +
# coord_flip() +
# labs(x = "rational for exclusion", y = "frequency") +
# scale_fill_brewer(palette = "Paired")
ggplot(evidence, aes(disturbance, fill = paradigm)) +
geom_bar(na.rm = TRUE) +
coord_flip() +
scale_fill_brewer(palette = "Paired") +
labs(y = "frequency")
#ggplot(evidence, aes(region, fill = paradigm)) +
# geom_bar(na.rm = TRUE) +
# coord_flip() +
# scale_fill_brewer(palette = "Paired") +
# labs(y = "frequency")
#ggplot(evidence, aes(data, fill = paradigm)) +
# geom_bar(na.rm = TRUE) +
#coord_flip() +
#scale_fill_brewer(palette = "Paired") +
#labs(y = "frequency")
#ggplot(evidence, aes(system, fill = paradigm)) +
# geom_bar(na.rm = TRUE) +
# coord_flip() +
# scale_fill_brewer(palette = "Paired") +
# labs(y = "frequency")
#ggplot(evidence, aes(goal, fill = paradigm)) +
#geom_bar(na.rm = TRUE) +
#coord_flip() +
#scale_fill_brewer(palette = "Paired") +
#labs(x = "outcome", y = "frequency")
#step 1 models####
#paradigm
derived.evidence <- evidence %>%
group_by(technique, data, region, disturbance, goal, paradigm) %>% summarise(n = n())
#active-passive split
m <- glm(n~paradigm, family = poisson, derived.evidence)
anova(m, test="Chisq")
## Analysis of Deviance Table
##
## Model: poisson, link: log
##
## Response: n
##
## Terms added sequentially (first to last)
##
##
## Df Deviance Resid. Df Resid. Dev Pr(>Chi)
## NULL 167 9.9147
## paradigm 1 0.045115 166 9.8696 0.8318
#region
m1 <- glm(n~paradigm*region, family = poisson, derived.evidence)
#m1
#summary(m1)
anova(m1, test="Chisq")
## Analysis of Deviance Table
##
## Model: poisson, link: log
##
## Response: n
##
## Terms added sequentially (first to last)
##
##
## Df Deviance Resid. Df Resid. Dev Pr(>Chi)
## NULL 167 9.9147
## paradigm 1 0.045115 166 9.8696 0.8318
## region 6 0.301367 160 9.5682 0.9995
## paradigm:region 6 0.213627 154 9.3546 0.9998
#outcome
m2 <- glm(n~paradigm*goal, family = poisson, derived.evidence)
#m1
#summary(m1)
anova(m2, test="Chisq")
## Analysis of Deviance Table
##
## Model: poisson, link: log
##
## Response: n
##
## Terms added sequentially (first to last)
##
##
## Df Deviance Resid. Df Resid. Dev Pr(>Chi)
## NULL 167 9.9147
## paradigm 1 0.045115 166 9.8696 0.8318
## goal 6 0.240941 160 9.6287 0.9997
## paradigm:goal 4 0.301480 156 9.3272 0.9897
#even split between active and passive evidence by all key categories
A summary of sort process using PRISMA.
library(PRISMAstatement)
prisma(found = 1504,
found_other = 5,
no_dupes = 1039,
screened = 1039,
screen_exclusions = 861,
full_text = 178,
full_text_exclusions = 101,
qualitative = 77,
quantitative = 66,
width = 800, height = 800)
Check data and calculate necessary measures.
#all data
#data <- read_csv("data/data.csv")
#data <- data %>%
#mutate(lrr = log(mean.t/mean.c), rii = ((mean.t-mean.c)/(mean.t + mean.c)), var.es = ((sd.t^2/n.t*mean.t^2) + (sd.c^2/n.c*mean.c^2)))
#other effect size estimates
#library(compute.es)
#data <- data %>%
#mutate(d=mes(mean.t, mean.c, sd.t, sd.c, n.t, n.c, level = 95, #cer = 0.2, dig = 2, , id = ID, data = data))
#check metafor
#data from ag and grazing studies that examined restoration in drylands
data <- read_csv("data/data.csv")
mydata<- data %>%
filter(disturbance == c("agriculture","grazing")) %>%
filter(!notes %in% "couldnt extract data")
mydata <- mydata %>%
mutate(lrr = log(mean.t/mean.c), rii = ((mean.t-mean.c)/(mean.t + mean.c)), var.es = ((sd.t^2/n.t*mean.t^2) + (sd.c^2/n.c*mean.c^2))) #use only lrr now
Explore summary level data of all data. Explore aggregation levels that support the most reasonable data structure and minimize non-independence issues.
#evidence map####
require(maps)
world<-map_data("world")
map<-ggplot() + geom_polygon(data=world, fill="gray50", aes(x=long, y=lat, group=group))
#map + geom_point(data=paperdata, aes(x=long, y=lat)) +
#labs(x = "longitude", y = "latitude") #render a literal map, i.e. evidence map, of where we study the niche in deserts globally
#add in levels and color code points on map####
#map of 178 articles
map + geom_point(data=data, aes(x=long, y=lat, color = paradigm)) +
scale_color_brewer(palette = "Paired") +
labs(x = "longitude", y = "latitude", color = "")
#map of selected articles, agriculture and grazing disturbances
map + geom_point(data=mydata, aes(x=long, y=lat, color = paradigm)) +
scale_color_brewer(palette = "Paired") +
labs(x = "longitude", y = "latitude", color = "")
#aggregation####
se <- function(x){
sd(x)/sqrt(length(x))
}
data.simple <- mydata %>%
group_by(study.ID, paradigm, technique, measure.success) %>%
summarise(n = n(), mean.lrr = mean(lrr), mean.rii = mean(rii), mean.var = mean(var.es))
main.data <- mydata %>%
group_by(study.ID, paradigm, intervention, outcome) %>%
summarise(n = n(), mean.lrr = mean(lrr), mean.rii = mean(rii), mean.var = mean(var.es))
#EDA data####
simple.data <- mydata %>% group_by(study.ID, paradigm, technique, measure.success) %>% summarise(mean.rii = mean(rii), error = se(rii))
simple.data <- na.omit(simple.data)
parad.data <- mydata %>% group_by(study.ID, paradigm) %>% summarise(mean.rii = mean(rii), error = se(rii))
parad.data <- na.omit(parad.data)
tech.data <- mydata %>% group_by(study.ID, technique) %>% summarise(mean.rii = mean(rii), error = se(rii))
tech.data <- na.omit(tech.data)
success.data <- mydata %>% group_by(study.ID, measure.success) %>% summarise(mean.rii = mean(rii), error = se(rii))
success.data <- na.omit(success.data)
#active
active <- mydata %>%
filter(paradigm == "active")
#viz for aggregation####
disturbance.data <- mydata %>%
group_by(study.ID,disturbance, paradigm) %>%
count()
disturbance.data2 <- disturbance.data %>%
group_by(disturbance,paradigm) %>%
count()
ggplot(na.omit(disturbance.data2), aes(disturbance,nn, fill=paradigm))+
geom_bar(stat = "identity") +
coord_flip(ylim=0:44) +
scale_fill_brewer(palette = "Blues")
intervention.data <- active %>%
group_by(study.ID,intervention, paradigm) %>%
count()
intervention.data2 <- intervention.data %>%
group_by(intervention,paradigm) %>%
count()
#ggplot(na.omit(intervention.data2), aes(intervention,nn, fill=paradigm)) +
#geom_bar(stat = "identity") +
#coord_flip(ylim=0:44) +
#scale_fill_brewer(palette = "Blues")
technique.data <- mydata %>%
group_by(study.ID,technique, paradigm) %>%
count()
technique.data2 <- technique.data %>%
group_by(technique,paradigm) %>%
count()
technique.data3 <- technique.data2 %>%
group_by(technique) %>%
count()
ggplot(na.omit(data.simple), aes(technique, n, fill = paradigm)) +
geom_bar(stat = "identity") +
coord_flip() +
scale_fill_brewer(palette = "Paired")
ggplot(na.omit(data.simple), aes(measure.success, n, fill = paradigm)) +
geom_bar(stat = "identity") +
coord_flip() +
scale_fill_brewer(palette = "Paired")
#better
ggplot(main.data, aes(intervention, n, fill = paradigm)) +
geom_bar(stat = "identity") +
coord_flip() +
scale_fill_brewer(palette = "Paired") +
labs(fill = "")
ggplot(main.data, aes(outcome, n, fill = paradigm)) +
geom_bar(stat = "identity") +
coord_flip() +
scale_fill_brewer(palette = "Paired") +
labs(fill = "")
Exploratory data analyses to understand data and QA/QC using Rii.
Meta and conventional statistical models to explore relative efficacy.
Step 5. Synthesis stats
#p-value meta
#library(metap)
#all data for metas but cleaned for na's
#all data for meta cleaned
mdata.all <- mydata %>%
filter(!is.na(lrr)) %>%
filter(!is.na(var.es)) %>%
filter(!is.na(n.t)) %>%
filter(!is.na(p)) %>%
filter(!is.na(intervention)) %>%
filter(is.finite(lrr))
#active only data for meta cleaned up
mdata <- mydata %>%
filter(paradigm == "active") %>%
filter(!is.na(lrr)) %>%
filter(!is.na(var.es)) %>%
filter(!is.na(n.t)) %>%
filter(!is.na(p)) %>%
filter(!is.na(intervention)) %>%
filter(is.finite(lrr))
#aggregated data for metas var estimated with central tendency
#note - could also bootstrap mean variance here instead of arithematic mean
simple.mdata <- mdata %>%
group_by(intervention) %>%
summarise(lrr = mean(lrr), var.es = mean(var.es), n = mean(n.t))
simple.mdata.2 <- mdata %>%
group_by(intervention, outcome) %>%
summarise(lrr = mean(lrr), var.es = mean(var.es), n = mean(n.t))
simple.mdata.var <- mdata %>%
group_by(intervention) %>%
summarise(lrr = mean(lrr), var.es = se(var.es), n = mean(n.t))
simple.mdata2.var <- mdata %>%
group_by(intervention, outcome) %>%
summarise(lrr = mean(lrr), var.es = se(var.es), n = mean(n.t))
#metas with p-values####
#schweder(mdata$p)
#sumz(p, data = mdata)
#mdata %>%
#split(.$intervention) %>%
#purrr::map(~sumz(.$p))
#sumlog(mdata$p)
#metas with meta package on effect size measures####
library(meta)
#all data non-aggregated
#m <- metagen(lrr, var.es, studlab = ID, byvar = intervention, data = mdata)
#summary(m)
#funnel(m)
#metabias(m)
#forest(m, sortvar = intervention)
#t-tests if different from 0
tmu <- function(x){t.test(x, mu = 0, paired = FALSE, var.equal=FALSE, conf.level = 0.95)}
mdata.all %>%
split(.$paradigm) %>%
purrr::map(~tmu(.$lrr)) #note this uses arithmetic means not estimated means from random effect models
## $active
##
## One Sample t-test
##
## data: x
## t = 3.4074, df = 359, p-value = 0.0007302
## alternative hypothesis: true mean is not equal to 0
## 95 percent confidence interval:
## 0.1316402 0.4909925
## sample estimates:
## mean of x
## 0.3113163
##
##
## $passive
##
## One Sample t-test
##
## data: x
## t = -5.8407, df = 171, p-value = 2.559e-08
## alternative hypothesis: true mean is not equal to 0
## 95 percent confidence interval:
## -0.5034850 -0.2491308
## sample estimates:
## mean of x
## -0.3763079
#Propose we only use random effect models because of the diversity of studies
m <- metagen(lrr, var.es, studlab = ID, byvar = paradigm, data = mdata.all)
summary(m) #active is positive and passive is negative so should not mix
## Number of studies combined: k = 496
##
## 95%-CI z p-value
## Fixed effect model -0.1348 [-0.1348; -0.1348] -111701369.65 0
## Random effects model -0.0139 [-0.0855; 0.0578] -0.38 0.7045
##
## Quantifying heterogeneity:
## tau^2 = 0.2648; H = 12688836.57 [12688835.55; 12688837.59]; I^2 = 100.0% [100.0%; 100.0%]
##
## Test of heterogeneity:
## Q d.f. p-value
## 79698253852777280.00 495 0
##
## Results for subgroups (fixed effect model):
## k 95%-CI Q
## paradigm = active 333 -0.1348 [-0.1348; -0.1348] 79607575146295024.00
## paradigm = passive 163 -0.1574 [-0.1574; -0.1574] 90524299681441.44
## tau^2 I^2
## paradigm = active 0.2649 100.0%
## paradigm = passive 0.3666 100.0%
##
## Test for subgroup differences (fixed effect model):
## Q d.f. p-value
## Between groups 154406800798.70 1 0
## Within groups 79698099445976464.00 494 0
##
## Results for subgroups (random effects model):
## k 95%-CI Q
## paradigm = active 333 0.1947 [ 0.1033; 0.2862] 79607575146295024.00
## paradigm = passive 163 -0.3537 [-0.4882; -0.2192] 90524299681441.44
## tau^2 I^2
## paradigm = active 0.2649 100.0%
## paradigm = passive 0.3666 100.0%
##
## Test for subgroup differences (random effects model):
## Q d.f. p-value
## Between groups 43.68 1 < 0.0001
##
## Details on meta-analytical method:
## - Inverse variance method
## - DerSimonian-Laird estimator for tau^2
#metabias(m)
#should do a t-test right now of paradigm is diff from 0
#interventions
m <- metagen(lrr, var.es, studlab = ID, byvar = intervention, subset = paradigm == "active", data = mdata.all)
summary(m)
## Number of studies combined: k = 333
##
## 95%-CI z p-value
## Fixed effect model -0.1348 [-0.1348; -0.1348] -111668536.95 0
## Random effects model 0.1947 [ 0.1033; 0.2862] 4.17 < 0.0001
##
## Quantifying heterogeneity:
## tau^2 = 0.2649; H = 15484891.12 [15484889.86; 15484892.37]; I^2 = 100.0% [100.0%; 100.0%]
##
## Test of heterogeneity:
## Q d.f. p-value
## 79607575146295024.00 332 0
##
## Results for subgroups (fixed effect model):
## k 95%-CI
## intervention = vegetation 221 -0.1348 [-0.1348; -0.1348]
## intervention = soil 75 0.7272 [ 0.7269; 0.7274]
## intervention = water addition 37 0.0000 [-0.0014; 0.0014]
## Q tau^2 I^2
## intervention = vegetation 79607575088050608.00 0.2649 100.0%
## intervention = soil 2099026.81 0.0509 100.0%
## intervention = water addition 2.34 0 0.0%
##
## Test for subgroup differences (fixed effect model):
## Q d.f. p-value
## Between groups 56145406.47 2 0
## Within groups 79607575090149632.00 330 0
##
## Results for subgroups (random effects model):
## k 95%-CI
## intervention = vegetation 221 0.1403 [ 0.0378; 0.2428]
## intervention = soil 75 0.4617 [ 0.3633; 0.5601]
## intervention = water addition 37 0.0000 [-0.0014; 0.0014]
## Q tau^2 I^2
## intervention = vegetation 79607575088050608.00 0.2649 100.0%
## intervention = soil 2099026.81 0.0509 100.0%
## intervention = water addition 2.34 0 0.0%
##
## Test for subgroup differences (random effects model):
## Q d.f. p-value
## Between groups 91.80 2 < 0.0001
##
## Details on meta-analytical method:
## - Inverse variance method
## - DerSimonian-Laird estimator for tau^2
#metabias(m)
m <- metagen(lrr, var.es, studlab = ID, byvar = intervention, subset = paradigm == "passive", data = mdata.all)
summary(m)
## Number of studies combined: k = 163
##
## 95%-CI z p-value
## Fixed effect model -0.1574 [-0.1574; -0.1574] -2736465.45 0
## Random effects model -0.3537 [-0.4882; -0.2192] -5.15 < 0.0001
##
## Quantifying heterogeneity:
## tau^2 = 0.3666; H = 747523.89 [747522.42; 747525.37]; I^2 = 100.0% [100.0%; 100.0%]
##
## Test of heterogeneity:
## Q d.f. p-value
## 90524299681441.44 162 0
##
## Results for subgroups (fixed effect model):
## k 95%-CI
## intervention = vegetation 54 0.2697 [ 0.2697; 0.2697]
## intervention = grazing exclusion 12 0.4327 [ 0.4327; 0.4328]
## intervention = soil 97 -0.7666 [-0.7666; -0.7666]
## Q tau^2 I^2
## intervention = vegetation 775623742262.70 0.0075 100.0%
## intervention = grazing exclusion 6686.32 <0.0001 99.8%
## intervention = soil 11085228014467.45 0.1148 100.0%
##
## Test for subgroup differences (fixed effect model):
## Q d.f. p-value
## Between groups 78663447918024.95 2 0
## Within groups 11860851763416.47 160 0
##
## Results for subgroups (random effects model):
## k 95%-CI
## intervention = vegetation 54 0.0434 [-0.0088; 0.0956]
## intervention = grazing exclusion 12 0.3481 [ 0.3412; 0.3550]
## intervention = soil 97 -0.4835 [-0.5740; -0.3929]
## Q tau^2 I^2
## intervention = vegetation 775623742262.70 0.0075 100.0%
## intervention = grazing exclusion 6686.32 <0.0001 99.8%
## intervention = soil 11085228014467.45 0.1148 100.0%
##
## Test for subgroup differences (random effects model):
## Q d.f. p-value
## Between groups 446.83 2 < 0.0001
##
## Details on meta-analytical method:
## - Inverse variance method
## - DerSimonian-Laird estimator for tau^2
#metabias(m)
#use random effects mean and var estimate to plot
#do t-tests here too
#outcomes
m <- metagen(lrr, var.es, studlab = ID, byvar = outcome, subset = paradigm == "active", data = mdata.all)
summary(m)
## Number of studies combined: k = 333
##
## 95%-CI z p-value
## Fixed effect model -0.1348 [-0.1348; -0.1348] -111668536.95 0
## Random effects model 0.1947 [ 0.1033; 0.2862] 4.17 < 0.0001
##
## Quantifying heterogeneity:
## tau^2 = 0.2649; H = 15484891.12 [15484889.86; 15484892.37]; I^2 = 100.0% [100.0%; 100.0%]
##
## Test of heterogeneity:
## Q d.f. p-value
## 79607575146295024.00 332 0
##
## Results for subgroups (fixed effect model):
## k 95%-CI Q
## outcome = soil 143 -0.1348 [-0.1348; -0.1348] 79607530173242944.00
## outcome = plants 108 -0.1990 [-0.1990; -0.1989] 68146638.29
## outcome = animals 11 1.0512 [ 0.9693; 1.1331] 0.05
## outcome = habitat 71 0.5577 [ 0.5577; 0.5577] 3395032696.46
## tau^2 I^2
## outcome = soil 0.2649 100.0%
## outcome = plants 0.3814 100.0%
## outcome = animals 0 0.0%
## outcome = habitat 0.0739 100.0%
##
## Test for subgroup differences (fixed effect model):
## Q d.f. p-value
## Between groups 41509872733.84 3 0
## Within groups 79607533636422272.00 329 0
##
## Results for subgroups (random effects model):
## k 95%-CI Q
## outcome = soil 143 0.1906 [ 0.0791; 0.3020] 79607530173242944.00
## outcome = plants 108 0.5161 [ 0.2783; 0.7538] 68146638.29
## outcome = animals 11 1.0512 [ 0.9693; 1.1331] 0.05
## outcome = habitat 71 -0.4266 [-0.5776; -0.2755] 3395032696.46
## tau^2 I^2
## outcome = soil 0.2649 100.0%
## outcome = plants 0.3814 100.0%
## outcome = animals 0 0.0%
## outcome = habitat 0.0739 100.0%
##
## Test for subgroup differences (random effects model):
## Q d.f. p-value
## Between groups 344.41 3 < 0.0001
##
## Details on meta-analytical method:
## - Inverse variance method
## - DerSimonian-Laird estimator for tau^2
#metabias(m)
m <- metagen(lrr, var.es, studlab = ID, byvar = outcome, subset = paradigm == "passive", data = mdata.all)
summary(m)
## Number of studies combined: k = 163
##
## 95%-CI z p-value
## Fixed effect model -0.1574 [-0.1574; -0.1574] -2736465.45 0
## Random effects model -0.3537 [-0.4882; -0.2192] -5.15 < 0.0001
##
## Quantifying heterogeneity:
## tau^2 = 0.3666; H = 747523.89 [747522.42; 747525.37]; I^2 = 100.0% [100.0%; 100.0%]
##
## Test of heterogeneity:
## Q d.f. p-value
## 90524299681441.44 162 0
##
## Results for subgroups (fixed effect model):
## k 95%-CI Q tau^2
## outcome = habitat 45 0.2697 [ 0.2697; 0.2697] 775623742237.89 0.0075
## outcome = plants 21 0.4327 [ 0.4327; 0.4328] 6703.49 <0.0001
## outcome = soil 97 -0.7666 [-0.7666; -0.7666] 11085228014467.45 0.1148
## I^2
## outcome = habitat 100.0%
## outcome = plants 99.7%
## outcome = soil 100.0%
##
## Test for subgroup differences (fixed effect model):
## Q d.f. p-value
## Between groups 78663447918032.59 2 0
## Within groups 11860851763408.83 160 0
##
## Results for subgroups (random effects model):
## k 95%-CI Q tau^2
## outcome = habitat 45 0.0091 [-0.0455; 0.0637] 775623742237.89 0.0075
## outcome = plants 21 0.3483 [ 0.3415; 0.3552] 6703.49 <0.0001
## outcome = soil 97 -0.4835 [-0.5740; -0.3929] 11085228014467.45 0.1148
## I^2
## outcome = habitat 100.0%
## outcome = plants 99.7%
## outcome = soil 100.0%
##
## Test for subgroup differences (random effects model):
## Q d.f. p-value
## Between groups 464.58 2 < 0.0001
##
## Details on meta-analytical method:
## - Inverse variance method
## - DerSimonian-Laird estimator for tau^2
#metabias(m)
#t-tests for outcomes diff from 0 but you can see using the 95% CI what is different
#check metafor for interactions?? in these big models or are we ok??
#brainstorm on how to explore?? techniques
#save but cut all this.
#m.study <- metagen(lrr, var.es, studlab = study.ID, byvar = intervention, data = mdata)
#summary(m.study)
#funnel(m)
#metabias(m)
#forest(m, sortvar = intervention)
#aggregated data
#m1 <- metagen(lrr, var.es, studlab = intervention, data = simple.mdata)
#summary(m1)
#funnel(m1)
#metabias(m1)
#forest(m, sortvar = intervention)
#different var estimate
#m2 <- metagen(lrr, var.es, studlab = intervention, data = simple.mdata.var)
#summary(m2)
#funnel(m2)
#metabias(m2)
#forest(m, sortvar = intervention)
#m3 <- metagen(lrr, var.es, studlab = intervention, byvar = outcome, data = simple.mdata.2)
#summary(m3)
#funnel(m)
#radial(m3)
#metabias(m2)
#forest(m, sortvar = intervention)
#m4 <- metagen(lrr, var.es, studlab = intervention, byvar = outcome, data = simple.mdata2.var)
#summary(m4)
#funnel(m)
#radial(m4)
#metabias(m2)
#forest(m, sortvar = intervention)
#forest(m1, sortvar = intervention, fontsize = 10, xlim = c(-2, 2))
#forest(m2, sortvar = intervention, fontsize = 10, xlim = c(-2, 2))
#forest(m3, sortvar = intervention, fontsize = 10, xlim = c(-2, 2))
#forest(m4, sortvar = intervention, fontsize = 10, xlim = c(-2, 2))
#t-tests for lrr diff from 0
#mdata %>%
#split(.$intervention) %>%
#purrr::map(~sumz(.$lrr))
#effect sizes plots
#need better ci estimates
#ggplot(simple.mdata, aes(intervention, lrr)) +
# ylim(c(-2,2)) +
# geom_point(position = position_dodge(width = 0.5)) +
# labs(x = "", y = "lrr") +
# coord_flip() +
# geom_errorbar(aes(ymin=lrr-var.es, ymax=lrr+var.es), width=.05, position = position_dodge(width = 0.5)) +
# geom_hline(yintercept = 0, colour="grey", linetype = "longdash")
#ggplot(simple.mdata.2, aes(intervention, lrr, color = outcome)) +
# ylim(c(-2,2)) +
#geom_point(position = position_dodge(width = 0.5)) +
# labs(x = "", y = "lrr", color = "") +
# coord_flip() +
#geom_errorbar(aes(ymin=lrr-var.es, ymax=lrr+var.es), width=.05, position = #position_dodge(width = 0.5)) +
#geom_hline(yintercept = 0, colour="grey", linetype = "longdash")
#ggplot(mdata, aes(lrr, color = intervention)) +
#geom_freqpoly(binwidth = 0.5, size = 2) +
#xlim(c(-10, 10)) +
#labs(x = "lrr", y = "frequency", color = "") +
#geom_vline(xintercept = 0, colour="grey", linetype = "longdash") +
#scale_color_brewer()
#ggplot(mdata, aes(lrr, fill = intervention)) +
#geom_dotplot(binwidth = 1) +
# xlim(c(-10, 10)) +
# labs(x = "lrr", y = "frequency", fill = "") +
# geom_vline(xintercept = 0, colour="grey", linetype = "longdash") +
# scale_fill_brewer()